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Abstract 

We report evidence of a superconducting instability (of T lg symmetry) in the infinite- 
ly Anderson lattice in the presence of crystal fields of cubic symmetry. We assume 
a lattice of 4/ sites, each with a total angular momentum of J = 5/2 that is split 
by crystal fields into a low-lying doublet of symmetry and an excited quartet of 
Tg symmetry. Slave Bosons on the 4/ sites create and destroy 4/° configurations 
and Lagrange multipliers at each 4/ site enforce the occupancy constraint due to 
the infinite Coulomb repulsion. Quasiparticle interactions are due to exchange of 4/ 
density fluctuations, which are represented by fluctuations in the slave Bosons and 
Lagrange multipliers. We use the so-called analytic tetrahedron method to calculate 
the dressed (to order 1/N) Boson Green functions. In weak couping, the exchange 



of the dressed Bosons gives rise to a superconducting instability of Ti g , xy[x 2 — y 2 ), 
symmetry. The A lg , "s-wave", channel has strongly repulsive interactions and hence 
no pairing instability. The T 2g channel exhibits weakly repulsive interactions. Average 
quasiparticle interactions in the E g , x 2 — y 2 , 3z 2 — r 2 , channel fluctuate strongly as a 
function of the number of tetrahedra used to calculate the Bosonic Green functions, 
lending only weak evidence for an instability of E g symmetry. 

PACS No. 74.70.Tx 
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I. INTRODUCTION 



This paper is concerned with the effects of crystal electric fields on quasiparticle in- 
teractions in the Ce based heavy Fermion superconductor CeCu2Si2 (T c ^0.6 K|T[). 
In general, the heavy Fermion materials are examples of systems exhibiting strong 
correlations among the constituent particles 0] ||, [|J, [[J, |J, 0. Generally, the 
compounds are comprised of intermetallics and rare-earth or actinide atoms (such 
as uranium or cerium) with a strong on-site Coulomb repulsion. This large electro- 
static energy arises from the extremely localized nature of the 4f or 5f wavefunctions 
in the solid and markedly influences the electron occupation at these "rare-earth" 
sites. When hybridization between a rare-earth electron and a conduction electron 
is allowed, the physics of this strong interaction is communicated to the solid at 
large, giving rise to a metal of strongly correlated, interacting electrons. In such a 
system one might expect to find a ground state manifesting collective properties of 
the coupled rare-earth and conduction electrons, e. g. superconductivity or mag- 
netism. Indeed such ground states are seen. There are also heavy Fermion systems 
that apparently retain a metallic state down to zero temperature. To-date there are 
six known heavy Fermion superconductors, all containing either cerium or uranium. 

Our work is based on the infinite-U Anderson lattice, the details of which we shall 
discuss later. For the experts, we mention here that we use slave Boson operators to 
create or destroy 4/° configurations on the Ce sites, thereby avoiding the cumbersome 
Hubbard operators in the hybridization piece of the Hamiltonian f|,[f|,[nj. There is 



also a Lagrange multiplier to enforce unit occupancy of the 4/ multiplets at each Ce 
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site. Our work is novel in that we also include, at the Ce sites, crystal electric fields 
of cubic symmetry, which has the effect of splitting the spin-orbit coupled (J = 5/2) 
multiplet into a doublet (of T 7 symmetry) and a quartet (of r 8 symmetry). We take 
the T 7 doublet to be the ground multiplet, with a crystal field splitting, A C ef, to the 
T 8 quartet that is much larger than the Kondo temperature of the low- lying doublet, 
T o7 . {Acef > T o7 ) 

Previous theoretical work has focused on understanding the heavy Fermion com- 
pounds mainly through the SU (N) version of the periodic Anderson model |TT| , [0 , |L3| , 
• In the SU(N) model, the 4/ multiplet is iV-fold degenerate, and the (plane- 
wave) conduction bands are assumed iV-fold degenerate as well. The matrix element, 
V(k), for hybridization between a conduction electron and a 4/ electron, is taken to 
be isotropic in momentum space. 

Within the SU(N) model, Lavagna, Millis, and Lee [fL5|J , Auerbach and Levin jJEJ, 
and Houghton, Read, and Won [ 16] have studied quasiparticle interactions due to 
the exchange of 4/ density fluctuations. The lowest order diagrams contributing to 
the interactions are of order 1/N, where N is the 4/ multiplet degeneracy. In Ce, 
in the absence of crystal field splitting, the low- lying J = 5/2 multiplet is six- fold 
degenerate (N = 6). So it seems reasonable to truncate the diagrams at order 1/N. 
Lavagna, Millis, and Lee found such a spinless density exchange yielded a ci-wave 
superconducting instability in the spin-singlet pairing channel. 



F. C. Zhang and T. K. Lee[T7] have performed a more realistic calculation (as far as 
heavy Fermion compounds are concerned) by including spin-orbit coupling at the Ce 
sites and by returning to the two-fold degenerate conduction states. They included an 
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anisotropic hybridization matrix element, of the Coqblin-Schrieffer form||18||, between 
conduction and 4/ electrons . These spin-orbit coupled ions are assumed to sit in an 
overall spherically symmetric "host" (as in a jellium model). Unlike in the SU(N) 
model, in the even-parity pairing channel, Zhang and Lee found no superconducting 
instabilities of s-wave, d-wave, or g-wave symmetry. It is worth noting here that the 
mean field quasiparticle energy bands Zhang and Lee found are the same as those 



calculated by Zou and Anderson who used the KKR scheme and included spin- 
orbit coupling on the 4/ sites. 

We have looked at density-fluctuation induced quasiparticle interactions in the 
presence of crystal electric field splitting of the Ce J = 5/2 multiplets in the lat- 
tice. We find the anisotropy due to cubic symmetry qualitatively and quantitatively 
alters the interactions in comparison to the results of Zhang and Lee. In fact, we 
find evidence for a superconducting instability of T l9 (xy(x 2 — y 2 )) symmetry in the 
even-parity pairing channel. We also find weaker evidence for an E g (x 2 — y 2 ) pair- 
ing instability. We find no instability in the A\ g (the "s-wave" of cubic symmetry) 
channel, which is not surprising, given the strong Hubbard U repulsion built into the 
Anderson model. 

The details of our calculation are presented here as follows. In Sec. II, we briefly 
discuss experimental evidence for the existence of crystal electric fields of cubic sym- 
metry in CeCu2Si2 . In Sec. Ill, we introduce the Hamiltonian for the infinite- 
U Anderson lattice in the presence of crystal fields. The Hamiltonian formalism 
shall be retained throughout this paper, as opposed to functional integral techniques 



12|,|14],[15|,[20|. In Sec. IV, we report the mean field properties of our Hamil- 
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tonian, including the quasiparticle energies and states. In Sec. VI, we discuss our 
calculation of the dressed slave Boson Green functions, which include the effects of 
particle-hole excitations in the hybridization coupled conduction and 4/ electron sys- 
tem. We also explain (Sees. VII and VIII) our use of the so-called analytic tetrahedron 
method for performing the complicated three dimensional Brillouin zone integrals that 
arise in the evaluation of the particle-hole diagrams. Finally, in Sees. X and XI, we 
present our results for the superconducting instabilities in the presence of cubic sym- 
metry, which are based on the Fermi surface average of the quasiparticle-quasiparticle 
scattering amplitude. 

II. CRYSTAL ELECTRIC FIELDS 



The first experimental evidence of crystal electric fields in CeCu2Si2 came from 
the inelastic neutron scattering data of Horn et al. [fJl] . The data show a clear peak at 
an energy of 31.5 meV (~ 360K) and a weaker peak at approximately 12 meV. Since 
there was also some structure at 12 meV in the reference material LaCu2Si2 (which 
has no 4/ electrons), it seems likely that the lower energy peak is due to phonons. 
The data could be fit, however, by assuming a crystal field structure in which the 
6-fold degenerate J = 5/2 multiplets are split into three doublets, with excitation 
energies given by the measured energies of the inelastic peaks. Such a three doublet 
structure is exactly what one would expect for a tetragonal crystal, like CeCu 2 Si 2 . 

Later, measurements of the specific heat at high temperatures |[22| , however, did 



not agree with this interpretation of the multiplet structure. Bredl et al. plotted the 
4/ contribution to the specific heat of CeCu 2 .2Si2 by subtracting off the corresponding 
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data for LaCu2.2Si2- They found they could not fit the results with two Schottky 
peaks, as would be expected for three crystal field doublets. Instead, the excited 
magnetic states behaved as if there were a four-fold degenerate multiplet about 360 
K above a ground state doublet. This idea is also nicely corroborated by Bredl 
and co-workers' plot of the entropy as a function of temperature. The rather quick 
rise of the entropy from a value of i?ln2 to i?ln6 at about 300 K would suggest a 
quartet structure rather than two doublets separated in energy by about 100K. Such 
a structure would be appropriate in the presence of cubic symmetry. 

Further evidence for effective cubic symmetry in CeCu2Si2 comes from the dc 
susceptibility measurements by Steglich and co-workers on single crystal CeCu2Si2 
samples [p3| . They find a very weak anisotropy in the temperature dependence of 
Xdc for a magnetic field applied parallel or perpendicular to the c-axis of the unit 
cell. Furthermore, the isotropy of the slope of the upper critical magnetic field at T c , 
^(^c)lPli suggests that the effective mass of the quasiparticles is isotropic. Since 
the effective mass is technically a second-rank tensor, it can be isotropic only in 
the presence of cubic symmetry. In light of all this information, it is reasonable to 
assume that the small peak in the neutron data at 12 meV is not due to crystal field 
excitations and that there is effectively cubic symmetry at the cerium sites. 

In the presence of cubic symmetry, the J=5/2 multiplet is split into a doublet of 
T 7 symmetry and a r 8 quartet |p5 |. The crystal field-split states \Ta), where a labels 



the degenerate states for a given multiplet, are a linear combination of the eigenstates 
of the z-component of the total angular momentum, |m), where —5/2 < m < 5/2: 

|r 7 ,+l) = -wf|-5/2) + w||3/2) (1) 



|IV,-l> = -^|5/2> + y |l - 3/2> (2) 
|r 8 ,+2) = y||5/2) + y||-3/2) (3) 

|r 8 ,i) = |i/2> (4) 
|r 8 ,-i) = |-i/2) (5) 

|r 8 , -2} = y|| - 5/2} + y||3/2). (6) 

From the neutron scattering and specific heat data we know that the ground 
multiplet is the doublet. Thus, without addressing the true microscopic source of the 
crystal fields, we deduce a multiplet structure as shown in figure 1, where Aqef labels 
the size of the crystal field splitting. In general, equations p] - |6] can be expressed as 

|T, a) = ^2c ram \m), (7) 

m 

where the coefficients cr am can be read directly from the equations. 

We can use equation |7| to construct the matrix element, Vraa(k), for hybridization 
between a crystal field state with quantum numbers T, a and a plane wave conduction 
state with crystal momentum k and spin a. In the case of a single ion in the full 
J=5/2 manifold, one may use the form derived by Coqblin and Schrieffer 



-> / 4-ir / 7 — 2mn 

VUk) = y Y^ 3(xVoh] J 14 — Y i™~%( k ), (8) 

(-5/2 < m < 5/2), er = ±1 is the (pseudo)spin index, and the angular dependence is 
in the spherical harmonic. V Q k denotes the dependence of the hybridization strength 
on the magnitude of the momentum, which will be important only near the zone 
center (i.e. near |fc|=0). In fact, we can write V Q k = V g(k), where g(k) is a function 



of I A; I that goes to zero at the zone center like k 3 . V Q is the bare hybridization strength. 
Because the crystal field split 4/ states are just linear combinations of the J = 5/2 
states, the hybridization matrix elements in cubic symmetry are 

5/2 

Vr aa {k) = c r am V mcT {k). (9) 

m=-5/2 

Figure 2 represents the hybridization process pictorially. 

III. HAMILTONIAN 

We start with a few words on the important energy scales of CeCu2Si2 , based 



on the discussion of Kang et al.\Wj on the electron spectroscopic data available as 
of 1990. They analyzed their own data of the Ce 3d x-ray photoelectron spectrum 
(XPS) and 4/ bremsstrahlung isochromat spectrum (BIS); they also analyzed Ce 
4/ resonant photoelectron data (RESPES) of Parks et al. ||27|| . Kang and co-workers 
calculated the appropriate one-electron spectra from the impurity Anderson model, 
which showed reasonably good agreement with the data. On the basis of such a 
calculation, the authors claim that Coulomb energy for double occupation of a Ce 4/ 
site is U « 7eV. 

The large Coulomb energy for CeCu2Si2 prompts us to take the limit in which U 
goes to infinity, thereby forbidding hybridization processes that give rise to 4/ 1 — > 
4/ 2 valence fluctuations. This is a technical simplification for us, but even though 
U is indeed very large by solid-state physics standards, it may be that the physics 
of finite U is crucial to the understanding of heavy Fermions. Recently Dan Cox|[28|| 
has proposed the Quadrupolar Kondo Model (or Two-Channel Kondo Model) as an 
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explanation of the superconductivity and of the possible non-Fermi liquid behavior in 
uranium based heavy Fermions. The physics needed to get non-Fermi liquid behavior 
derives from the existence of enough independent channels of conduction electrons to 
overcompensate the effective 4/ "spin". A finite-?/ version of the Anderson model, 
as applied to a Ce impurity, can exhibit a two-channel Kondo effect. The crystal 
field split multiplet structure here is the key. Taking the Tj doublet as the low-lying 
multiplet for the 4/ 1 configuration, and a T3 (non-magnetic) doublet as the low lying 
multiplet for the 4/ 2 configuration (as is appropriate for cubic symmetry), a two- 
channel Kondo effect is possible. Conduction states of Tg symmetry can mix the 4/ 1 
and 4/ 2 configurations f29[| , and conduction states of Tj symmetry can mix the 4/ 1 



and 4/° configurations. If the hybridization of Tg conduction states is stronger than 
hybridization of T 7 conduction states, the physics of the two-channel Kondo effect will 
determine the low temperature properties of the model. It may be that the finite-U 
Anderson model is a better starting point for heavy fermion systems. Nevertheless, 
the infinite U limit is a reasonable simplification (at least for CeCu2Si2 ) of an already 
difficult problem and warrants study in its own right. 

Our Hamiltonian can be written as a combination of terms: H = H c + Hf + 
H mix + H constraint . The kinetic energy of the conduction electrons is given by, 

ku 

where 



h 2 k 2 

2m 



— "A**, C 11 ) 



is the plane-wave dispersion. The zero of energy for this calculation will be taken 
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with respect to fi Q , the chemical potential of the conduction electrons in the absence 
of hybridization. 

The 4/ electron site energy is 



Hf=J2 *fLjSra> ( 12 ) 



'RTa 

RTa 



where R is the site index in real-space, and Ta are the crystal field quantum numbers. 
The operator f^ Fa destroys a 4/ 1 configuration at lattice site R, in which the crystal 



field state Ta is initially occupied. From the photoemission data| 26| , we take the 
energy of the doublet to be -2.0 eV, i.e. E7 =-2.0 eV. From the inelastic neutron 
scattering data^TJ and the high temperature specific heat|22|, we take the Tg level 
to lie 360 meV above the T 7 level, 

E 8 = E 7 + Acef = -1.964 eV. (13) 

The hybridization, or mixing, term is 

(14) 



s kaRTa 



The operator is a slave Boson creation operator, which creates a 4/° configuration, 
or hole, at lattice site R. N s is the number of lattice sites. This combination of con- 
duction, 4/, and slave Boson operators was first applied to the Anderson model by 
Barnes||, and was later reintroduced by Coleman |10|]- H mix contains only Bosonic 



or Fermionic operators, and so Wick's theorem is applicable. Use of the more cum- 
bersome Hubbard projection operators, 

X 0Va = \0){Ta\, (15) 
11 



would make Feynman diagrammatic procedures invalid, since the Hubbard operators 
do not obey standard commutation relations. The cubic symmetry is reflected in the 
structure of the anisotropic function Vr QCT (/c), which has a significantly different k 
dependence than the Coqblin-Schrieffer form V mcT (k). This new anisotropy will affect 
quasiparticle interactions differently than in the case of spherical symmetry. 

Finally, there is the constraint term, which is introduced with a Lagrange multi- 
plier, insuring that the total occupancy (Fermions plus Bosons) at the Ce sites 
is unity, 

Hconstraint = ]T l\ & [f^J R Ta + blb & " 0) ■ (16) 
R 

Note that Q=l is the physically meaningful value in this case. 

There are a few technical points we need to mention here. As discussed by Read 
and Newns || in the functional integral approach to the infinite-U Anderson impurity, 
it is useful to write the slave Boson operator as the product of a modulus and a phase 
factor, which in the lattice problem picks up a site index, R, 

h = (17) 
A gauge transformation for the / operators, 

then absorbs all the phase factors. In the functional integral approach, furthermore, 
the Lagrangian corresponding to our Hamiltonian (equations [T2|, [TJ], and |TB|) 
contains an (imaginary) time derivative of the slave Bosons, 

r[J) Or ' 
12 



which, from equation [T7|, introduces a factor of the phase velocity, 9g, into the La- 
grangian. Traditionally, this phase velocity is absorbed into the Lagrange multiplier, 
elevating it to the status of a dynamical field, z'A^(t). In this paper, since we have 
used the Hamiltonian formalism, we remark that, in the limit of static slave Bosons, 
the Lagrange multiplier plays the same role as it does in the functional integral 
formalism We have a 2x2 matrix Green function for the slave Boson and La- 
grange multiplier, which will have elements composed of averages over the following 
combinations of fields: ss, sX, and AA. 

Our calculation is based on previous 1/N calculations of the type applied to the 
SU(N) model, where N is the degeneracy of the 4/ multiplet. In the SU(N) case, 
in order to have a well defined Kondo temperature in the limit of large N, it was 
necessary to assume the bare hybridization strength, V a , scaled like l/\/jV[II]. We 
do the same here by defining a rescaled hybridization matrix element, 



V raa (k) = ^N r V raa (k), (18) 

where we assume that Vr aa (k) is of order 1. Note that N r is the degeneracy of the T 
crystal field multiplet. Also, because the hybridization matrix element and the Boson 
operator, s^, appear together in H mix , we define a scaled Boson operator, 



where we assume that s^ r is also of order 1. 



IV. MEAN FIELD APPROXIMATION 



13 



In this section, we discuss the properties of our Hamiltonian at mean field level, 
where we assume both the Bose operator and the Lagrange multiplier are uniform in 
space. In this limit, the Hamiltonian in k-space takes the following form: 

Hmf = E £*4a + E e ^lvM 



ka kTa 



E 



SorV ra(7 (k)ciJ l : ra + H.c 



koTa 

N, 



+ -fY l N v i\ [fov-Q^)- (20) 
where er = Er + i\ is the shifted energy of the T multiplet. We assume that q Q r is 
of order 1, but technically, when it comes down to getting numerical results, we know 
that q Q r = 1/iVr. Note that s Q r and A Q are the mean field values of the Bose operator 
and Lagrange multiplier, respectively (s Q r = s / y/ N?) ■ 

It is straightforward to diagonalize Hmf and obtain the quasiparticle states and 
energies. Because the r§ states are four-fold degenerate, there is one non-hybridizing 
quasiparticle band of T 8 symmetry. The secular equation which gives the quasiparticle 
energies E n ^ {n— 1,2,3 is the band index) is 

''o7l v 7aa\~ 6 o8 1 v 8aa | 



e 8 - E 



nk) 



5 2 IT/ I 2 Z 2 IT/ I 2 

it _ rp \ \ "* ^o7l v 7aa\ _ \ -* VR, 



0, (21) 



( e 7~ E nk) as ( 6 8-^„fc) 

where the sum J2a 7 (s) means only the states in the r 7 ( 8 ) doublet are summed. Al- 
though technically it is possible to find analytic solutions of equation |2T], the analytic 
expressions are too cumbersome to be useful; thus we calculated the roots numerically. 
The three quasiparticles bands are plotted in figure 3 along two different directions 
in the cubic Brillouin zone. Figure 3 shows that along the axes of the zone, for ex- 
ample along the TX direction, the states of T-j symmetry can not hybridize with the 
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conduction states. This means that the matrix elements Vj acr (k) vanish along these 
special directions and that there is no gap between the first and second quasiparticle 
bands (see figure 3(b)). Such behavior has been discussed by Mart in [J50]. 

The quasiparticle states can, quite generally, be written as a combination of plane 
wave and crystal field states, 

. r . s _t7* (t\ i 

(22) 



\QsJ = Mk) 



\ ka ) -2s p — l ra ) 



where \ka) is a plane wave state. The anisotropic normalization function is 



A 2 n (k) = 



§ 2 r\Vr a *(k)\' 



(er - E. 



-1 



(23) 



For given crystal field quantum numbers, Fa, it is possible to sum over the pseudo- 
spin variable a in equation 23. For convenience, we define a function 



(k) = J2 V raAk)V a rw(k). 



(24) 



Then the sum over the pseudospin variable in equation ^3| is just the diagonal element 
[iTaTa{k)- All possible (non-zero) forms for the function //r<*r' are shown in Table 
1, where the dependence on the magnitude of k has been divided out. Using the results 
of this table, it is possible to arrive at the following expression for the normalization 
function, which is valid at any point in the Brillouin zone: 



A 



T o7 
SoVo 



I _|_ 2 / To7_ 
3 "t" 3 \ T o8 



2y^ 

9 



1 



T 8 



Y AQ {k) + M[Y^{k)+Y^ 4 {k){ 



(25) 



Figure 4 shows a plot of equation ^5] along the equator of a spherical Fermi surface, 
with denoting the azimuthal angle measured with respect to a coordinate axis. The 
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extremely sharp variations near the axes represent the vanishing of the hybridiza- 
tion matrix elements. We call the six points where the Brillouin zone axes intersect 
the Fermi surface "hot-spots"; they must be handled with care when averaging the 
quasiparticle scattering amplitude over the Fermi surface. 

We would like to remark that, in this model, the existence of "hot-spots" is a man- 
ifestation of the lowering of the symmetry below spherical. In spherical symmetry, 
there is a sum rule for the matrix elements that renders the normalization function, 
A 2 , isotropic in k-space[|U|],|T7]]. In our case, the matrix elements can not be sim- 



plified to an isotropic function. This means that "hot-spots" can occur anywhere in 
the Brillouin zone where some subset of the hybridization matrix elements, Vr arT (k), 
vanishes. This is a rather general statement that relies only on the symmetry being 
lower than spherical and should not be unique to cubic lattices. 

We can relate the fermionic creation and destruction operators in the original basis 
(the c and / operators) to creation and destruction operators in the quasiparticle 
basis, Qt na , Q ln(J . We find that 

^ = J2Mk)Q^ (26) 

n 

ncr fc r -C n £ 

These expressions are useful in constructing the two-quasiparticle scattering ampli- 
tude. 

In the last part of this section, we discuss the self-consistency equations that 
arise when one demands that the free energy in the mean field approximation be an 
extremum with respect to the Bose fields s Q and i\ Q . Since we have diagonalized the 
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mean field Hamiltonian, equation |20], we can write it in terms of the quasiparticle 
energies and operators: 

\ 

kncr ' 



Hmf = E E n% Q\ nrr Q %na + ^A - g or ) . (2? 



The corresponding mean field free energy has the form 

Fmf = y E ^A (s 2 or - qor) ~ j E ln(l + e"^), (29) 

where /3 is the inverse temperature. Requiring that OFmf I 'di\ =0 and OFmf /ds Q =0 
yields the equations 

1 1 r9P 

5 IWfc-*r) + ^ 2/(^)^ = 0, (30) 

and 

1 f) F 

2iA „ s „ + - £/(^)^ = 0. (31) 

fccrra 

In equations ^ and |3l| f(E n %) is the Fermi function evaluated at the quasiparticle 
energy B r. We also need an equation to fix the chemical potential of the quasipar- 
ticles, n, which depends on the total number of electrons (conduction electrons, n c , 
and f electrons, rif) per unit cell 

ntotai = n c + n f = — E f{ E n %)- ( 32 ) 

s kan 

These three coupled integral equations, when solved self-consistently, give the 
shifted 4f multiplet energies, ep, the value of the Bose field, s , and the quasiparticle 
chemical potential, \l. The input parameters are n tota i, the quasiparticle filling factor, 
the bare hybridization strength, V Q , and the conduction electron filling factor, n c . The 
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zero of energy is always measured relative the the chemical potential of the conduction 
electrons, fi Q . For n tota i =2, the lowest quasiparticle band is completely filled, and 
the system is a Kondo insulator. We have consistently used n tota i =1.5, which insures 
that we have a metal. 

To get numerical self-consistent solutions, we found it necessary to write these 
equations in terms of energy integrals with the appropriate density of states, 

ka 

where iV(£) is the density of states per spin for the unhybridized conduction elec- 
trons and dQ is an element of solid angle. Note that in the SU(N) model described 
previously, the spin degeneracy would contribute a prefactor of iV (instead of 2) in 
equation For free electrons in three dimensions the density of states is proportional 
to the square root of £. 

To proceed, we make one approximation. We assume that surfaces of constant 
energy for the quasiparticle states are spherically symmetric. Near the zone center, 
this is exactly correct, and there is no approximation at all. Near the zone boundary, 
the equal-energy surfaces become distorted from spheres due to the constraints of T 7 



symmetry. It is important to note, however, that in equations |3(J, |3l], and [32], the 
strongest angular dependence comes from the anisotropic matrix elements, which we 
treat exactly. That is, we believe the angular dependence of the quasiparticle bands 
is not as important as that of the hybridization matrix elements. For example, near 
an axis of the Brillouin zone, the T 7 matrix elements are going to zero. So, even if the 
quasiparticle energies surfaces are distorted from spheres near the axes, the sensitivity 
of the self-consistency equations to this distortion would be lessened by the presence of 
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the small V? a o- terms in the numerator. Thus, it should be a reasonable approximation 
to treat the mixing matrix elements as having all the angular dependence. 

Using the secular equation for the quasiparticle band energies, equation El], we 
can calculate all the necessary derivatives of the quasiparticle energies found in equa- 



tions |30] and After averaging the anisotropic matrix elements over the Fermi 
surface, we are left with the following three equations to be solved self-consistently: 

slvl r » .„ „ r i 5 



EMs 2 r-Qov)+2 S -^ I" dE N(£(E)) 

P 3 J-D 

2i\ Q - 2^ f dE N(£(E)) 

3 J-D 



+ 



(e 7 - Ef (e 8 - E)' 
5 



+ 



ntotai-2 I JE N(£(E)) 



e 7 -E e 8 -E 
5 

+ 77 



0, and 



d b o v ok 



0, (34) 

(35) 
(36) 



3(e 7 -£) 2 S^-E) 2 . 
where — D is the energy at the bottom of the lowest quasiparticle band. 

We have taken the limit of zero temperature to arrive at equations ^4|-0. As 
a consequence, only the first quasiparticle band, E\, contributes at mean field; for 
simplicity, we have dropped the band subscript "1". We consider two different sets of 
solutions, corresponding to the conduction electron filling factors of n c =0.5 (which 
we call set (a)) and n c =0.8 (which we call set (b)). We define a Kondo temperature 
in the lattice for both T 7 and T 8 multiplets by 



T r = e r - fJ,. 



(37) 



The reader is reminded that /i is the quasiparticle chemical potential. The motivation 
for defining the Kondo temperature as the difference between the shifted multiplet 
energy (er) and the quasiparticle chemical potential, is that in the SU(N) model this 
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difference has exactly the same structure as the Kondo temperature for the impurity 
problem. That is, in the SU(N) model, one finds that PI] 



e- fi = De-^ NN ^ v °, (38) 

where D is the half bandwidth, N(0) is the (assumed flat) conduction electron density 
of states, Ef is the unshifted 4/ multiplet energy, and V is the bare hybridization 
strength. In table 2, we see that both parameter sets (a) and (b) have approxi- 
mate Kondo temperatures (for the Tj doublet) of 10 K, which, based on the neutron 



scattering quasielastic linewidth[^T|], is a reasonable estimate for CeCu2Si2 . Note, 
also, that the Kondo temperature for the T 8 quartet is dominated by the crystal field 
splitting, 

T o8 = T o7 + Acef « 370 K, (39) 

where T j/T s ^0.027. Figure 5 plots the lowest quasiparticle bands for both mean 
field parameter sets (a) and (b). The top of the first quasiparticle band is just below 
the shifted T 7 energy, e-j. The quasiparticle chemical potential cuts through the flat 
part of the first band in both cases, giving rise to a very large quasiparticle density 
of states at the Fermi surface. 

V. EFFECTIVE MAGNETIC MOMENT 

In this section we discuss our results for the effective magnetic moment of the 
quasiparticle states when averaged over the Fermi surface. This calculation was mo- 
tivated by the work of Zou and Anderson ||19|| , who wished to explain how the small 



Wilson ratios, R, seen in heavy Fermion superconductors (for CeCu 2 Si 2 , R ~0.5||) 

20 



could be reconciled with the presence of strong spin fluctuations. If spin fluctuations 
are indeed the source of the superconductivity in heavy fermion compounds Q 



, one 



would expect a Wilson ratio greater than unity. Zou and Anderson have claimed that, 
even in the presence of strong spin fluctuations, the Wilson ratio could be reduced due 
to the anisotropic hybridization between the conduction and the / electrons. They 
calculated the quasiparticle states, including spin-orbit coupling, from the relativistic 
KKR equation. Their states are identical to those found via a mean-field approxi- 



mation to the infinite-U Anderson lattice[ll|, [|T^], [33], or via Gutzwiller projection 
techniques |34] . 



In spherical symmetry, for a total angular momentum of J = 5/2, the quasiparticle 

states are 



| for) 



\Qkna) = M\k\) , 

m=-5/2 

where the normalization for the lowest energy band is 



5/2 s V (k) 

E°o v ma v. / i \ 
T- ^l m ) 



E n k) 



(40) 



A\{\k\) 



1 + 



i 2 V 2 

o r o 



'lk 



(41) 



In the preceding equations, E nk is the quasiparticle energy for band n and e is the 
energy of the shifted 4/ states. Note that the normalization function is isotropic and 
depends only on the magnitude of the wavevector. For the case of J = 5/2 there 
is a four-fold degenerate band of energy e, since two of the six total / states have 
hybridized with the conduction states. 

Zou and Anderson calculated the magnetic moment of the quasiparticle states in 
the lowest band, which we denote as 



(42) 



21 



where the magnetic moment operator is 

fi z = l z + 2s z = gjJ z , (43) 

gj being the g-factor for the total angular momentum, J. The effective magnetic 
moment is the average of \i a z a ' over the spherical Fermi surface, 



dk 
An 



(44) 



and the result is /igjy=1.16/x^, where //g is the Bohr-magneton. Note that the bare 
moment for a f 1 state is /i=2.54/i B , so that hybridization has reduced the size of the 
moment. 

The Pauli susceptibility for the quasiparticles is a function of the effective magnetic 
moment and the renormalized density of states at the Fermi surface, N(0) (iV(0) = 
—N(0), where m* is the quasiparticle effective mass), 

XPauii = 2/4 ff N(0). (45) 

Substituting this value for the Pauli susceptibility into the expression for the Wilson 
ratio gives 

Rreduced = (^) = ^Ja' ( 46 ) 

In this way, it would be possible to have a large Stoner factor, (1 + F£ )~ 1 , as befits 
the presence of strong spin fluctuations, but the reduced magnetic moment could still 
account for a small experimentally observed Wilson ratio. 

Cox[^5j, Zhang and Lee[^], and Aeppli and Varma|37| pointed out the need to 



include the Van Vleck susceptibility in the calculation. The value for x i n equation |46 



must be the total susceptibility, as seen by experiment. In the case of the uniform 
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susceptibility, the Van Vleck contribution represents direct transitions from the Fermi 
energy in the lowest band to higher energy bands. The dominant contribution must 
come from transitions to the four-fold degenerate band lying at an energy T above 
the Fermi energy, where T is the Kondo temperature in spherical symmetry. The 
Van Vleck susceptibility then has the structure 

Xvv(?=0) = 2 ,l E f^ msJhMl , (47 ) 

where the quasiparticle density of states N(0) ~ 1/T D . The small energy denominator 
(for CeCu2Si2 , T Q w 1 meV) makes this interband contribution to the susceptibility 
comparable to the contribution from the Pauli susceptibility. Cox and Zhang and Lee 
showed that the total susceptibility was given by 

Xtotai = XPauli + Xvv = \g]J{J + l)/i|[2JV(0)], (48) 

which depends on the bare magnetic moment. Thus the inclusion of the Van Vleck 
susceptibility would alter Zou and Anderson's argument. 

We now wish to discuss what happens when crystal electric fields of cubic sym- 



metry are included. The quasiparticle states are given by equation |22], where the 
local crystal field states are orthonormal to each other. Note that in this case the 
dispersionless band lies at an energy Aqef above the Fermi energy, where the crys- 
tal field splitting is much larger than the Kondo temperature. In order to calculate 



the effective magnetic moment, we need the expectation value of equation |42] for the 



lowest quasiparticle band. We find that we can write equation [42] in the form 

nf{k)=A\{k) 



mgjCr' a 'mSor'V^ a , a (k)Vraa'(k)s rCTc 



F'aTam ~ E lk)(^ ~ E lk) 
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(49) 



Using equation [0], we find the average moment is 

vlubiceff = 0.583/4, ( 50 ) 
for mean-field parameter set (a). For a free Ti doublet the average moment is 

2 25 2 /r 1 \ 

^7 = 4gPB- ( 51 ) 

The effective magnetic moment in cubic symmetry is slightly larger than the value for 
a free moment. This was expected by CoxfSH], who predicted the effective moment 
would have the structure 

(52) 



2 2 



, m 

1 + a< 



where a is a prefactor that could be as big as about 10. The contributions of order 



m/m* come from the six regions on the Fermi surface which intersect the axes 
of the cubic Brillouin zone. These are the so-called "hot-spots", which we have 
already mentioned. At these points on the Fermi surface, the plane wave conduction 
states and the localized T-j states can not hybridize. Thus we expect the g-factor 
at these "spots" to go back to the free electron value of 2, with the result that the 
contribution to /i e // from the "hot-spots" increases the effective moment above that 
of a free moment of T 7 symmetry. For parameter set (a), we have 

™ = = 1.33 x 10- 

m* \s V J 

Equation is applicable for the quoted values of /j^^ and fij if a=9.23, which is a 
reasonable value. 

Finally, we substitute fi^jj- into the expression for the Wilson ratio, 

R =f^ym (53) 

V/W 7(0)' 1 ; 
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where x(0) is the low temperature susceptibility, 7(0) is the linear coefficient of specific 
heat, ks is Boltzmann's constant, and fiB is the Bohr magneton. Using x(0)=0.019 
emu/mole (for a magnetic field along the c axis of the tetragonal unit cell) |38|j , and 
7(0)=1000 mJ/mole-K 2 |l| gives R = 2.38. It is possible that experimentally quoted 
results of approximately 0.5 for the Wilson ratio || are too small since the bare 
effective moment of n 2 e jj =2.54/i^ was used in equation |52|. We also find it interesting 



to note that for the two-channel Kondo impurity, the Wilson ratio is approximately 
2.6. 

VI. HYBRIDIZATION DRESSED PROPAGATORS 



With the eigenstates of the mean field Hamiltonian written in the undiagonalized 
basis, there are three fermionic Green functions: Gr a r'a'(the f Green function); G^the 
conduction Green function); and Grw^he off-diagonal, or mixing, Green function). 
All three are defined below in terms of Fock space operators: 

GW&t) = -(Trh ra (r)fl, al (0)), (54) 

G a (k,r) = -(T T c %a (r)cl(0)), (55) 

G Faa (k,r) = -<T^r>)4 CT (°)>> ( 56 ) 
where T r is the imaginary time ordering operator. It is easy to see, diagrammatically, 
how these hybridization dressed Green functions can be calculated. Figures 6,7, and 
8 show the expansions for Gr a r a ', G a , and Graa, respectively. In terms of a complex 
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frequency, z, the mean-field Green functions are: 



z — er 



G a (k,z) = 

Graa(k,z) = 



_L (z - e 7 )(z - e 8 )s or VT aa (k)V r > a > a (k)s or i 

z ~ e r' (z - (z - E 2U ) (z - E^) 

(z - e 7 )(z - e 8 ) 



SoTV^ aa (k)(z-6 T *) 



(z - E$){z - E#){z - Eg)' 
These propagators will be used to calculate the dressed bosonic Green functions. 



(57) 
(58) 

(59) 



VII. FLUCTUATIONS 



In this section, we shall calculate the Green functions for the Bosonic fields in 
our Hamiltonian, namely s^ r and The Bosonic Green functions can be dressed 
(through the terms H mix and Hf) by particle-hole excitations of the hybridized con- 
duction and 4/ electron systems. To proceed, we write the slave Boson and Lagrange 
multiplier as follows: 

hr = 5 °r% + <% r , (60) 

i\ n = iX 6 nfi + i6X n , (61) 

where <5§£ r and SX^ represent fluctuations away from the (self-consistent) mean field 
values. 

It is easy to see how the particle-hole excitations dress the Bosonic Green functions 
by writing the full Hamiltonian in k-space, 



V ft 

"ka~ko 1 jL^t J kVa 
ka kk'Va 



E rkk> + iX k-k' 



fie' 



k'Ta 
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kk'Taa 



kl 



+ ^T, N ^ X k E r UkT S kr - Qr 



(62) 



Substituting from equations |60] and ^T], equation £2] can be written in two pieces, one 
representing the mean field approximation (which we have solved), and the second 
piece coming from the fluctuations in the Bosonic fields, Ss^ T and SX^. The bare Boson 
Green functions come from the terms in the constraint, the last line of equation ^2| 
which are quadratic in the fluctuating fields. Since this involves terms of the form 
Ss^ r SSj: r , 5Sj: r 5\j:, and SX^SXt, we write the bare propagator in a matrix form (in the 
static limit), 

D orr> = — 7T\ s n °rr', (63) 



where s G r = s Q / ^/N~r~. The specific elements of the matrix are D oss rr> = — ^i\ 8rr 
Dosxrv = -^fs or 5 r r>, and D oX \rr>=0- 

The dressed Boson Green function will then satisfy a matrix Dyson's equation, 



(64) 



The 2x2 self-energy matrix, tlrr'io), due to the particle- hole excitations, can be 
calculated by the usual Feynman diagrammatic techniques. 

To motivate the results for D^,, consider the diagonal component of the matrix 
Dyson's equation in Figure 9(a). Taking all contributions to the self-energy of order 
N-p would give a dressed propagator of order 1/Nr- In Figure 9(c), the x symbols at 
the corners of the diagrams represent scaled hybridization matrix elements, Viw(&)> 
which are of order one ((9(1)). Thus, since there is a sum over the degenerate states 
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a of the multiplet V and no factor of 1/iVr from the matrix elements to cancel it, 
we say that the Boson self-energy diagram is of order N-p. Technically, there is not a 
free sum over the degeneracy label a of a given multiplet, because the hybridization 
matrix elements, Vp aa (k), are dependent upon a. (That is, we do not have an explicit 
factor of Nr after summing over a.) In the SU(N) model, the isotropic hybridization 
matrix element, V Q , is truly independent of the degeneracy label, m. Thus the sum 
over m yields a factor of N. In our case, even though we do not have explicit factors 
of Ny-, we assume that the closed Fermionic bubbles are the important diagrams, in 
analogy with the SU(N) model. 

Figure 9 shows all the unique self-energy diagrams in terms of the three hybridization- 



dressed Fermionic Green functions, equations |54|-|56|. Evaluation of theses diagrams 



leads to the following results for the elements of the (inverse) dressed Bose propagator: 

D^ rr ,(q) = 2I ssrr ,(q), (65) 



—s r + — 
2 s oT . 



Srv +H s \rr'(q), (66) 



DxxrAQ) =-£-- \hxrr{q), (67) 
where the momentum dependent functions are given by 

N s ££ W ( E lP ~ E lk) ( E 2k> ~ E lk) ( E 3k> ~ E lk) ( E 2k ~ E lk) 



i^k',k+q Jr %',fc-qJ ^gg-j 



(er-^)(er-^)(%-^)' 
W(g1 = L J- f( E ik)^- E ik) 2 ^- E ik? 

" " N s m , ( E 2k ~ E lk)( E 3k ~ E lk)( E W ~ E lk)( E 2k' ~ E lk)( E 3k> - E lk) 
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{k)fir a r'a'(k ; ) 



lk>M<r +S k>,k- g l (69) 



" ' N ' ££ w ( E lk' ~ E lk)( E 2k> ~ E lk)( E 3k> ~ E lk)( E 2k ~ E lk) 

x (e 7 ~ ^) 2 (e 8 - E^f 

In equations |68| - [7(], _E n £ are the quasiparticle band energies. The quasiparticle 
chemical potential was chosen to lie in the lowest band, such that f(E n ^) is nonzero 
only for n — 1. (Recall this is zero temperature calculation.) The anisotropic function 
firar'a'{k) is defined in equation 123. The parameters x-p and yr are numbers which 



depend on the mean field parameter set used. See Table 3 for the values of these 
parameters. 

Note the presence of the principal value integrals over the Brillouin zone. We are 
forced to evaluate the integrals numerically, as discussed in the next section. This is 
the most labor-intensive part of the calculation. 

VIII. ANALYTIC TETRAHEDRON METHOD 



We wish to evaluate the three dimensional principal value integrals that have 
arisen in the calculation of the Bosonic self-energy, i.e. the functions I ss rr'(q), 
I s \rr' (<?) , and I\\rr' (<?*) of equations |68| - [70| . To do so, we use a procedure developed 
originally to calculate the real part of spectral functions of the formf39|j , [f40f , fflfl , [f42f , 



$'(£)= i?e$(£) = P^£-^^. (71) 
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$' looks like the real part of a dynamic susceptibility. Handled numerically, the 
principal value nature of the sum in equation |7T] makes it non-trivial. In fact, if we 
tried to do the sum by simply evaluating the function at many points throughout 
the zone, and then multiplying by a weighting factor, we would find an essentially 
infinite variance. That is, upon averaging the result over many different mesh sizes, 



the variance of the mean would be huge compared to the mean- value it self [|43 



The analytic tetrahedron method |4l|,|4^| evaluates expressions like that of equa- 



tion |7T] by breaking up the Brillouin zone into tetrahedra, where the band energies E ? 
need to be known only at the four corners of a given tetrahedron. Any band energy 
for k inside a tetrahedron is interpolated from the energies at the corners. Within this 
assumption of linear interpolation, and also assuming the effective matrix element, 
M(k), is constant inside the tetrahedron, Rath and Freeman and independently 



Lingard showed that the principal value integral over the tetrahedron could be 
performed analytically. 

We write our functions I a brr>{o) i n the following form: 

T (m _P^ fiMM^IliMl s (72) 

labvr [q) - -Tf — E - _ E - — °k>,n±? y u ) 

s kk' 1^' lk 

where the matrix elements are a complicated function of k and q which come from 
equations |6R KM and 170. That is, direct comparison of equation 173 with equations I 



70] gives the structure of the matrix elements M a bTT'k, k' . Note that ab can represent 
ss, sX, or AA. The Fermi function, of course, restricts us to only the volume below 
the Fermi surface. Since the surfaces of constant energy are planar, we need to 
know the possible unique ways a plane can cut through a tetrahedron. It turns out, 
fortunately, that there are only three such ways; and furthermore, in each of the three 
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cases, the two subdivided volumes of the intersected tetrahedron (the volume above 
the Fermi surface and the volume below the Fermi surface) are themselves either a 
single tetrahedron or a composite of three tetrahedra. We never deal with anything 
but tetrahedra. So it is just a matter of geometry to find out which tetrahedra are 



beneath the Fermi surface and hence contribute to the integral of equation |72 . 

To evaluate our self-energies, we need the analytic result for the following integral 
over a tetrahedron: 



tetra 



E-.r, — E,r' 



J lk' 

where k' = k ± q. For simplicity of notation we define the energy differences 



lkV 



(74) 



where the "1" signifies the band index, and i — 1, 2, 3, 4 labels the four corners of the 
tetrahedron. The result for the integral / depends only on the V, and the volume of 



the tetrahedron. Using the notation of Rath and Freeman||l|] the results is: 

Vi I ^ 2 



'V? 
3w te i ( -p- ln 



V 2 , 



Yi 



+ 



— In 



V 
Vi 



(75) 



where Vtet is the volume of the tetrahedron, and where 



(76) 
(77) 
(78) 

We have studied equation [75] carefully in the limits where some of the Vs are equal 
to each other or are equal to zero. It turns out that there are quite a few such cases, 





(Vi 


- V 4 )(Vi 


-V 3 )(Vi 


-v 2 


D 2 = 


(V 2 


-V 4 )(V 2 


-V 3 )(V 2 


- Vi 


D 3 = 


(V 


-V 4 )(V 3 


-V 2 )(V 3 


-Vi] 



but for the sake of space, we do not tabulate them here [[44 
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Note that these three dimensional integrals can have an entire surface of poles in- 
side the Brillouin zone that must be handled properly. In contrast, a one-dimensional 
principal value integral can be regularized numerically by basically subtracting off 



the divergence [45] , and such a procedure is facilitated by the relatively small number 



of poles throughout the domain of integration. In our three-dimensional integral, if 
n mesh i s the number of sub-cubes inside the cubic Brillouin zone, then there are of the 
order of n^ es/l poles, which clearly gets large as n mes h increases, and it is no longer 
possible to regularize the integral in a simple way. The large quantity of work done 
by (mostly electronic structure) physicists in this field of numerical k-space sums, is 
indicative of the degree of complexity inherent to these problems [pj,[|7 . 



We now show the results for the self-energy functions I ss (q), Is\(<f), and I\\(q), 
which have been summed over the crystal field multiplet indices, 

Ug)=E ^f 1 "® (79) 



rr' x 
and 

Mg) = EW(?l (81) 

rr' 

For a given mesh parameter, n mes h, the total number of tetrahedra in the Brillouin 
zone is 8xn^ esft . The functions I ss , I s \, and I\\ are plotted as a function of the mesh 
parameter, n mes h, i n Figures 10-12 for mean field parameter set (a). The results for 
parameter set (b) are similar. We chose the momentum q = 0.5z, measured in units 
of it /a, where a is the lattice spacing. (We take a =3.89A.) The matrix elements are 
assumed constant inside a given tetrahedron. Since they are complicated anisotropic 
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functions, however, the matrix elements can vary a great deal throughout the Brillouin 
zone, and thus they cause fluctuations in the value of the I a b(q) as a function of n mes h- 
This is clearly visible in figures 10-12. 

It is clear that if we need convergence of the Green functions to several decimal 
places, it would require a mesh with n mes h > 80. In background, a calculation for 
n mes h = 40 takes about one hour on a DEC 5100. Since the computing time goes as 
n mesh> ^ 1S clear that a run with n mes h = 80 requires a long run time. Attempts at 
vectorization of the code were hindered by a plethora of logic statements required for 
the subroutine that performed the integral in equation |73|. Actually, this numerical 
integration is parallel in nature. Each of the v? mesh subcubes could, in principle, be 
integrated independently of the others, and at the end the net result would be the 
sum of the results from each sub cube. Naively, this is the kind of problem a parallel 
machine should be able to handle well. 

It is also interesting to note that calculations of the susceptibility or dielectric 
constant based on electronic structure data, generally are not performed for n mes h 
greater than about 30 [48]. The complicated nature of our matrix elements have 
forced us to push the procedure to very large values (by anyone's standards) of n mes h. 
We shall discuss in the next section, how our conclusions on pairing instabilities in 
the infinite-U Anderson lattice take into account this slow convergence of the Bosonic 
self-energy. 

IX. QUASIPARTICLE SCATTERING AMPLITUDE 

In this section, we present our analysis of the scattering amplitude TQ P (k,k'), 
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which allows for quasiparticle interactions via the exchange of 4f density fluctuations. 
The exchanged Boson will be represented by the (dressed to order 1/N) Boson matrix 
Green function calculated in Sec. VIII. The diagrams for the scattering amplitude 
are presented in Figure 13, where the straight lines represent the incoming (and out- 
going) quasiparticles and the wavy line is an element of the (matrix) Boson Green 
function, Drr'(q)- The large black circles represent quasiparticle vertices, "ffr(k, k') 
and jmixrikjk'), which are calculated by writing the full Hamiltonian in the quasi- 
particle basis: 

Hqp = E kK{^A n 4)Q\ na Q %nl(T 

kann 1 

+ E Vr(l k')Qt \E T 5 rk , + Q Pn , a , 



' kna 

kk'Tann'aa' 



+ _ E [IrnixTik^Q^Q^^sl ^ + H.C. 
kkTaann' a' 



H constraint- 



Recall that Q kna destroys a quasiparticle of momentum k, pseudospin a, and band 
index n. The vertex functions 7/ and j m i X thus come from the unitary transformation 



that diagonalizes Hmf (see equations ^6] and |27|) 

7/r< ' ' = («r " ~ E^l ' m 

wE _ P) s _K(k)A A Z) loT v^)v T «Ai) (84) 

e r - E n>%> 

As shown below, we project TQ P (k,k') onto states of cubic symmetry, 
the so-called cubic harmonics, where 77 labels the irreducible representations of the 
octahedral group O n . The product is then averaged over the Fermi surface, 
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A superconducting instability of symmetry 77 is signaled by a negative value of the 
corresponding average, T v . See Table 4 for a list of the cubic harmonics used in this 
calculation. The character table for the octahedral group O is presented in Table 5, 
where the irreducible representations are listed: Ai, A 2 , E, T 1; and T 2 . The group 
Oh follows from the group O by including inversions. This means the representations 
pick up a subscript " g" or "u" depending on if the they are even or odd, respectively, 
under parity. Because of the complexity in calculating the dressed bose propagators 



we have assumed a spherical Fermi surface for the average in equation |85|. We do not 
feel this is a weakness of the calculation for reasons discussed elsewhere |f49| . 

In this paper, we shall discuss only even parity pairing states. This restriction 
is based on the following experimental evidence for CeCu2Si2 : the need for strong 
Pauli limiting to fit the low temperature upper critical field data, if C 2(0)[^3|]; the 
reduced spin susceptibility below T c as measured by the 63 Cu Knight shift pUL and 
the observed T 3 temperature dependence of the nuclear-spin relaxation rate below 
T c [50]. The strong Pauli limiting actually only argues against equal spin pairing 
states. As Ueda and Rice showed [fnjl , in the presence of spin-orbit coupling, Pauli 
limiting is possible for pairing states of T lu or T 2u symmetry. (Both of these odd- 
parity states would have a gap with point nodes as opposed to line nodes on the 
Fermi surface.) These experimental facts put together, however, might be considered 
reasonable evidence for even-parity pairing in CeCu2Si2 . 

We have found it useful to study the properties of the scattering amplitude in 
two steps. First, by setting the functions I ss rr'{(l), I s \rr>{(]), and I\\rr'{(f) to zero, 
we simplify the problem considerably to that of two quasiparticles scattering via 
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exchange of a momentum independent Boson. In real space, this corresponds to a 
local interaction between the quasiparticles. We find these local interactions, when 
averaged over the Fermi surface, are substantially different from those calculated 
within the jellium model by Zhang and T. K. LeejlTj . 

In the case of cubic symmetry, inclusion of the functions I ss tt'(q), Isxrvio), and 
Ixxrv (o) is the computationally intensive part of this calculation. If there are strong 
local repulsive interactions in the rj pairing channel, then the only way to get a pairing 
instability (i.e. T v < 0), is to have the functions I S srr'(<l), hxrvio), an d I\xrr'(o), 
which represent the effect of non-local interactions, overcome the repulsion. Zhang 
and Lee discovered that in spherical symmetry these q-dependent contributions are 
too weak to overcome the local repulsions in the s, d, and g-wave pairing states. We 
find that in cubic symmetry (with crystal-field splitting) attractive nonlocal interac- 
tions can overpower local repulsions in the T\ g pairing channel, thus giving evidence 
for a T\g pairing instability. 

When we evaluate the diagrams of Figure 13, we can write the quasiparticle scat- 
tering amplitude in the even-parity (pseudospin singlet) channel as 

D ss vr{k' -k) + D ssTT ,(-k' -%)- i — — [D s xr{k' - k) + D sX r{-k' — %)) 

e r - E x i> \ J 

SoT ' f D sXr (k' - k) + D sXr (-k' - k')^ 

D xx (k' -k)+D xx (-k' -k))}, (86) 



x 



e r' - E w 

SoFSoF' 



(e T - E l% ,){e T , - E l% .) 
where the symbols mean the following: 
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*Vo-ra{k) = y/ NrV a ra is the hybridization matrix element between a plane wave 
conduction state of momentum k and spin a and a crystal field state with quantum 
numbers T and a. N-p labels the degeneracy of the F multiplet. 
•s Q r = s a / \J Np can be thought of as the mean field hybridization renormalization 
coefficient. At mean field level, the bare hybridization is renormalized (due to the 
constraint of only allowing hopping onto an empty 4/ site) to the value s D rV^ra- 
•er is the self-consistent, shifted mean field energy of the T multiplet. 
•E^ is the quasiparticle energy for the lowest, or first, band as a function of momen- 
tum. 

•A\(k) is the quasiparticle normalization function defined in equation 23[ and which 



is very strongly peaked at the six points where the Brillouin zone axes intersect the 
Fermi surface. 

•The components of the Bosonic Green function are: 

D ssTr (q,T) = (5s_gr(r)(5sgr/(0)); 

D sXr (q,r) = (5~s^(r)5X q iO)}; 

D xx (q,r) = (5X^)6X^0)). 

Note that we are taking the static limit of these Green functions. It became clear to us 
that including the frequency dependence of the boson Green functions was impossible, 
given the difficult numerical integrals encountered even in the static limit. 

When evaluated on the Fermi surface (at zero temperature), the band energies 
and E^, are set equal to the quasiparticle chemical potential, /i. Thus a term 
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like er — E^ becomes 

e r - n = T oF , (87) 

which is the Kondo temperature of the T multiplet. For the case of CeCu2Si2 , the 
Kondo temperature of the T 7 doublet is approximately 10 K, while for the T 8 quartet, 

T o8 = T o7 + A CEF = 370K. (88) 



This explains how to treat all the energy denominators in equation ^6. 

As mentioned, it is instructive to consider the so-called local limit of equation KB 
in which the self-energy functions I ss (o), and I\\(q) are set to zero, 

1 _ A?(fc)A^^0s o rs o r^ r * aCT (^0^Tw(^^r*av(^^ra(^) 
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where 



r osA = W^V rS ~ or + ^), (90) 
r s yJ 



SoT 

r oAA = E| L - (9i) 

The combination of normalization functions, A\(k)A\(k'), and the product of the 
four hybridization matrix elements is due to the anisotropic vertices j m ixr(k, k') and 
7/r'(fc, k'). The contribution from the Bose Green function is that which remains 
inside the square brackets in equation |S9|. 



The local scattering amplitude of equation |89| is much simpler to deal with than 
the full expression of equation |36j. Roughly, the functions I ss {q), Is\(q), and I\\(q) can 
be thought of as renormalizing the local interactions. In the next section, we present 
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our results for the Fermi surface averaged scattering amplitude, r^, both with and 
without the Bosonic self energy. 

X. RESULTS-LOCAL LIMIT 



We now present our results for the Fermi surface averaged (local) quasiparticle 
interactions, Ti oca[ ri , where r\=A lg , E g , T l9 , or T 2g . We have not included pairing 
states in the A 2g representation (See Tables 4 and 5) for two reasons: (1) the lowest 
order spherical harmonic present in the cubic harmonic <&A 2g is ^6m[H], and the 
relatively rapid variation throughout the Brillouin zone of the pair wavefunction would 
correspond to a pairing state with a high kinetic energy and hence should be less 
accessible than the pairing states labeled by the other representations; (2) the A 2g 
state would also require a finer mesh for the Fermi surface average than the one 
we have used and hence would further increase the (already considerable) overall 
computing time. 

Table 6 gives the ratio Ti oca[jri /T o7 for both mean- field parameter sets (a) and (b). 
Both parameter sets have Kondo temperatures of about 10 K. Also presented in table 6 
is the contribution to the Fermi surface average from the so-called "hot-spots" . For the 
hot-spot contribution, the normalization functions A\(k) and A\{k') in equation 



were approximated with delta functions in k-space that sampled only the six points 
on the Fermi surface intersected by the Brillouin zone axes. The integrated weight of 
the delta functions was chosen to equal the area under the peaks, as shown in figure 
4. Such a calculation gives us a feeling for the importance of these special points 
where the normalization function, A 1 (k), is rapidly changing. This contribution is 
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labeled by Ti oca i^ot in Table 6. We see from Table 6 that 

• The Ai g pairing channel has (by-far) the largest (repulsive) local interaction, but 
it is not dominated by what happens at the "hot-spots" . 

• In the Eg channel, the strong anisotropy of the normalization function, gives rise 
to a weakly attractive local interaction. In fact, the attractive local interaction here 
is dominated by the contributions from the hot-spots. We note from Table 5 that a 
pairing state of E g symmetry transforms like x 2 — y 2 or 3z 2 — r 2 , which have maxima 
along the directions of the axes. The normalization is also strongly peaked along the 
axes. Thus the E g states are greatly affected by A\(k). 

• In the T\g and T2 9 channels, the local interaction is weak and repulsive. Neither 
symmetry channels can "see" what happens at the "hot-spots", because both cubic 
harmonics, $Tig and $T2g vanish identically at these points. 

In Table 7 we present for comparison, the results for the Fermi-surface averaged 
scattering amplitude from Zhang and Lee's|17] jellium model calculation, where Leg- 



endre polynomials (P v ) play the role of the cubic harmonics. States described by r]=0 
are s-wave; rj—2 corresponds to d-wave; and ??=4 is g-wave. Zhang and Lee find all 
non-zero interactions are repulsive and of about the same strength. As we have said, 
they also found the inclusion of non-local interactions was not sufficient to overcome 
the local repulsions. In octahedral symmetry, however, the very weak local inter- 
actions in the E g , Ti g , and T 2g pairing channels make a superconducing instability 
likely. 

XI. RESULTS-INCLUDING NON-LOCAL INTERACTIONS 
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This section discusses our results for the Fermi surface average of the full quasipar- 
ticle scattering amplitude, Yq P {k, k'), as given in equation |8(| The functions I ss rr'(o)i 
hxrr'io), and Ixxrr'io) (see equations ^S|-^), contain the physics of the screened (by 
density fluctuations of the coupled conduction-4/ electrons) slave Bosons. The mo- 
mentum dependence represents the contribution of non-local quasiparticles interac- 
tions in real space. 



From equations |86] and [85], we see that we must evaluate the self-energy functions 
for all unique combinations of k±k'. This we do separately from the actual averaging 
process. We store the required values for I ss rr'(<l), hxrvio), and I\xrr'(<l) in a look-up 
table. 

We discussed the general convergence properties of the functions I ss rr'(q), Is\rr'{q) 
, and I\\rr'(<l) in Sec. VIII, where we saw that, due to the complicated effective matrix 
elements in cubic symmetry (see equation [T^), they are very slow to converge as a 
function of the Brillouin zone mesh size. The slow convergence prompted us to try the 



following line of attack. We have averaged the scattering amplitude (equations [86D 
over the Fermi surface with a fixed averaging mesh. We vary, however, the mesh for 
calculating the Bosonic Green function, as characterized by the parameter n mesh , and 
study the average interactions as a function of n mes h. With this procedure, it became 
clear that n mes h « 50 is a practical limit of the mesh size. Creating the look-up table 
for the Bosonic self-energies would require well over a week of runtime on a DEC5100 
for anything bigger. Thus, it is important to ask if we can make any conclusions 
about possible pairing instabilities for n mesh < 50. 

The values of Y n /T o7 for r\ = A lg and r\ = T 2g are plotted as a function of n mesh 
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in Figure 14 for parameter set (a) and in Figure 16 for parameter set (b). Also on 
the plots, the local value for the average are marked by horizontal lines for each 
representation (except A lg ). The results show the following: 

• In the Aig channel, the interactions are clearly repulsive. The non-local interactions, 
however, are attractive, since the average local repulsion is of the size Ti oca i t Ai g /T o7 =3 .54, 
and the inclusion of the momentum dependent Bosonic self-energy gives Tai 9 /T o7 ss 
1.2. Thus the non-local contribution has reduced the local repulsion by about a factor 
of two. 

• In the T 2g channel, but for a glitch at n mes h=23, the interactions are repulsive for 
parameter set (a), with T T 2 g /T o7 ^0.20. For parameter set (b) (Figure 16), however, 
it is less clear if the average interaction is attractive or repulsive. Taking the results 
for both parameter sets together, we believe there is no pairing instability of T 2g 
symmetry. 

The values of T rj /T o7 for 77 = E g and rj = Ti g are plotted in Figure 15 for parameter 
set (a) and in Figure 17 for parameter set (b). 

•The results for T lg show attractive interactions, with an average value T T i g /T o7 =- 
0.212±0.025 for parameter set (a), and =-0.226±0.049 for parameter set (b). 

• The interactions in the E g channel are less well behaved, but are attractive for 
25< n mesh < 47 (for parameter set (a)). These numbers point to a possible super- 
conducting instability of E g symmetry, but the results have not converged enough 
for us to be sure. As a rough guide, using the numbers for 25< n mes h < 47 would 
give (T) Eg /T o7 ^-0.163±0.054. The values for parameter set (b), however, have a 
surprisingly large fluctuation at n mes h—31, which makes it very difficult to make any 
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conclusion for n mes h < 35. It is clear that the E g states are the most sensitive to 
the fluctuations due to the anisotropic matrix elements in the Bosonic self-energies. 
From these results, we see that the T l9 pairing state is the most likely candidate for 
a pairing instability in cubic symmetry. 

XII. DISCUSSION 



As we have said, n mes h ~50 is a practical limit on the size of the Bosonic self-energy 
mesh that we could run on a local workstation. The appearance of fluctuations in 
the averaged scattering amplitude (Figures 14-17), which limit our conclusions about 
possible pairing instabilities in cubic symmetry, is not shocking. After all, we know 
that there are fluctuations in the Bosonic self-energies for n mes h ~80 (See figures 11- 
13). As we discussed, such variations are due to the anisotropic matrix elements in 
the principal value integrals over the Brillouin zone. The zone is divided into a large 
number of tetrahedra (8x(n mes(l ) 3 to be exact), and the matrix elements are assumed 
constant inside a given tetrahedron. If the matrix elements are sharply peaked in 
some region of the Brillouin zone, then it is easy to see how fluctuations can occur. 
If the mesh is constructed so that the matrix elements are evaluated very near a 
peak, then only a small shift in the mesh is required before the matrix elements will 
be evaluated at a point far down on the sides of the peak. Thus a small change in 
mesh size could result in a large change in the evaluated matrix elements. A similar 
problem can arise in finite size lattice problems, where large variations in results can 
persist up to very large system sizes |53|]. The solution for that particular problem is 



an average over boundary conditions. Our problem, unfortunately, has no such cure, 
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and we must live with reasonable conclusions from the data we are able to gather. 
To reiterate, we have been able to draw the following conclusions: 

• The best candidate for a pairing instability is in a state of T lg (xy(x 2 +y 2 )) symmetry, 
with T Tlg /T o7 =-0.212±0.025. 

• The Eg channel (x 2 — y 2 , 3z 2 — r 2 ) also shows weak signs of an instability, with 
r Eg /T o7 =-0.163±0.054. 

• Quasiparticle interactions are strong and repulsize in the A ig ("s-wave") pairing 
channel. 

• A pairing instability of T 2g (xy, xz, yz) symmetry also appears highly unlikely. 
Using the classic weak-coupling equation for the superconducting transition temper- 
ature within a given representation, rj, we find (in the T\ g channel) 

T c (r] = T lg ) = 1.13T o7 e To7/rT1 « « 0.09 K, (92) 

which is smaller than the measured T c of CeCu 2 Si 2 . We did not expect, however, 
such a weak-coupling calculation to give a quantitatively accurate value for T c , which 
leads us to say a few words about strong versus weak coupling results. 

In the case of heavy Fermions, the effective Fermi temperature is of the order 
of the Kondo temperature, which is also the energy scale of importance for the su- 
perconducting glue. That is, on physical grounds there is no reason to believe that 
only a very thin energy shell (thin compared to the Kondo temperature, T o7 ) about 
the Fermi surface is of importance for superconductivity. Thus, a strong-coupling 
calculation, including the energy dependence of the scattering amplitude, should be 
performed. Given the complexity of the static problem, however, including the dy- 
namics of the slave Bosons in the presence of crystal fields is not feasible. This does 
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not necessarily imply that the pairing instabilities found in the static problem have 
no meaning. We shall discuss, generally, why this is so. 

It is conventional wisdom [p^], that near the Fermi surface the quasiparticle self- 
energy for heavy Fermions is strongly frequency dependent but only weakly dependent 
on the magnitude of the momentum, \k\. This can be understood intuitively as follows. 
The characteristic energy scale for the quasiparticles is the Kondo temperature, T o7 , 
which is about 10K. The degeneracy temperature for a typical metal is Tp « 10,000K. 
The characteristic momentum, however, is set by the Fermi wavevector, fcp, which 
for CeCu2Si2 is the size of a typical metal. Thus, broadly speaking, we expect the 
quasiparticle self-energy, E, to behave (near the Fermi surface) as 

-"Of-1, 
dui \T o7 J' 

d(g u {d)' 

where D is the bandwidth of the conduction band. Then, since 

as <9£ 

it seems reasonable to ignore the momentum dependence of the self-energy. This 
intuitive result was reinforced by Millis and Lee[ll[], who found (in the SU(N) model) 
that the momentum dependence of the imaginary part of the conduction electron self- 
energy (at order 1/N) is very weak (going as 1/iV 2 ), while the frequency dependence 
goes as the inverse of the Kondo scale. 

If one accepts the dominance of the frequency dependence in the quasiparticle self- 
energy, then it seems reasonable to assume that including such dependence in a strong 
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coupling calculation (ala McMillan [55]) would serve only to reduce the transition 
temperature, T c . We do not believe that the sign of the average scattering amplitude, 
T^, would be affected by such frequency dependence. Thus, our conclusions about 
which pairing channels, rj, show a superconducting instability should not be changed 
as the result of a strong coupling calculation. 

Please note that we are not saying that the reduction of the transition temperature 
due to the frequency dependence of the residual quasiparticle interactions would be 
identical in structure to the case of electron-phonon coupling |55j. We can only say 
now that we expect T c to be reduced; we can not give an estimate of how large the 
reduction would be. 

The superconducting instabilities themselves appear to be based heavily on the 
underlying symmetry of the problem, which does not care if one performs a strong 
or weak coupling calculation. The fact that our estimated (weak-coupling) transition 
temperature for the T\ g pairing instability is smaller than the measured value for 
CeCu2Si2 , is not a surprise. The major purpose of this calculation has not been to 
give a precise numerical recipe for calculating the T c of heavy Fermion systems. We 
wished to study the importance of local, or "multiplet" , physics upon quasiparticle 
interactions. Thus the a weak coupling calculation should be a reasonable starting 
point. 
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TABLE CAPTIONS 



Table 1: Functional forms for the angular dependence of fir a r'a' f° r & U possible 
combinations of the crystal field quantum numbers T,a,T',a'. Note, that by time- 
reversal symmetry, firar'a'=^r*a*r'*a'* , where the * denotes the time-reversed pair. 
For example, ^71,71 = ^7-1,7,-1- Note, also, that V D r(|A;|) represents the dependence 

— * 

of the hybridization strength on the radial component of k. Any combination of 
quantum numbers not present in the table or not the time-reversed pair of quantum 
numbers in the table, will vanish upon summing over the pseudo-spin indices. 

Table 2: Self-consistent mean field parameter sets, labeled as (a) and (b). At 
mean field there are three coupled integral equations, which are solved self-consistently. 
The input parameters are: the bare hybridization strength, V Q ; the lower edge of the 
conduction electron band in the absence of hybridization, -D; the total number of 
electrons, n to tai = n cond + n f, per unit cell; the bare, or unshifted, energy of the T 7 
(EV) and the r 8 (E 8 = E 7 + Acef) multiplets; and the fixed crystal field splitting, 
Acef = 360 K. The self-consistent parameters which solve the equations are: the 
hybridization renormalization coefficient, s Q , where the mean field renormalized hy- 
bridization is s V ; the shifted T 7 (e 7 ) and T 8 (e 8 = € 7 + Acef) multiplet energies; the 
quasiparticle chemical potential, /i; and the Kondo temperature T o7 = e 7 — /i of the 
T 7 doublet. In both parameters sets (a) and (b), the total number of particles was 
fixed at n tota i=1.5. And the input parameters were chosen to give approximate Kondo 
temperatures of 10 K, i.e. T o7 ^10 K. All energies are measured relative the chemical 
potential in the absence of hybridization. In both parameter sets, A C ef=360 K, and 
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the unshifted T 7 energy is E 7 =-2.0 eV. 

Table 3: The parameters xr and yr are for both mean field parameter sets 
described in Table 2. 

Table 4: Realizations of the cubic harmonics, as linear combinations of the 
spherical harmonics Yi m . For each representation, 77, of cubic symmetry, the expansion 
was cut off after the lowest set of spherical harmonics with I > 0. 

Table 5: Irreducible representations of the octahedral group, O. To get the group 
O/j, we add inversions to the allowed symmetry operations, the result of which is that 
all representations pick up a subscript, g (for even parity) or u (for odd parity). 

Table 6: The local ("hard-core") quasiparticle scattering amplitude in the pres- 
ence of crystal electric fields, Tqp, averaged over a spherical Fermi surface. Results 
for both mean field parameter sets (a) and (b) are given. See table 5.3 for a discussion 
of the mean field parameters themselves. The first column gives the representations, 
labeled by 77, of the group Oh- The second and third columns are the averaged local 
scattering amplitudes (divided by the T 7 Kondo temperature, T o7 ) for the parameter 
sets (a) and (b), respectively. The fourth column lists (for parameter set (a)) the 
contribution to the average from the so-called hot-spots, where the Brillouin zone 
axes intersect the Fermi surface. In column four, in the E g channel, the attractive in- 
teraction seems to be due to these hot-spots; but in the A lg channel the hot-spots do 
not dominate, since the full Fermi surface average is large and positive. The average 
in the E g should be most sensitive to the hot-spots, since that is where the E g cubic 
harmonics have their maximum value. In the T ig and T 2g channels there is rigorously 
zero contribution from the hot-spots because the T ig and T 2g cubic harmonics vanish 
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at those points. 

Table 7: The local ("hard-core") quasiparticle scattering amplitude in the jellium 
model, as studied by Zhang and Lee[l]. When divided by the Kondo temperature 
and averaged over the spherical Fermi surface, these results are universal; there are 
no other parameters involved. In spherical symmetry, the pairing states are labeled 
by their relative angular momentum, /, with 1=0 corresponding the s-wave; 1=2, d- 
wave; and 1=4, g-wave. It is meaningful to compare Zhang and Lee's results with 
ours, because in the limit of spherical symmetry, where all the Bosonic propagators 
in equation 5.37 are replaced with their values in spherical symmetry, and where the 
normalization functions, A\(k) and A\(k'), are replaced by their isotropic values in 
spherical symmetry, our expression for Tqp (equation 5.73) gives exactly the same 
local interactions as Zhang and Lee. Note that all local interactions are repulsive in 
spherical symmetry. 
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FIGURE CAPTIONS 



Figure 1: Splitting of J = 5/2 multiplet into a T 7 doublet and a T 8 quartet due 
to crystal fields of cubic symmetry. Neutron scattering gives a splitting of about 360 
K. 

Figure 2: A schematic representation of the hybridization process, in which a 

— * 

conduction electron in the state \ka> jumps into an empty 4/ orbital of total angular 
momentum J = 5/2, and z-component m described by |m>. The matrix element is 
<m\V\kcr>= V ma (k). If we expand the conduction state in partial waves, then only 
the state with total angular momentum J = 5/2 and z-component m would be able 
to hybridize. Note that the energy of the J=7/2 excited multiplet is physically too 
large; it is shown thus merely to make the figure more readable. The actual size of 
the spin orbit coupling in CeCu 2 Si 2 , should be about an order of magnitude larger 
than the crystal field splitting. Note, also, that the energy splittings of the T 6 , T 7 , 
and T 8 states, which come from the crystal field splitting of the J=7/2 multiplet, 
are not accurate. The states are shown merely to inform the reader how the J=7/2 
multiplet decomposes in cubic symmetry. See, for example, Appendix D in reference 
[40]. 

Figure 3: Schematic of the quasiparticle bandstructure, showing the shifted crys- 
tal field multiplet energies {e^ and eg), the chemical potential (//), and the lower band 
edge (-D). (a) the k values range from the zone center (the Y point) to the intersection 
of the k x axis with the cubic Brillouin zone boundary, (b) From the zone center to 
the intersection of a cube diagonal with the Brillouin zone boundary. 

55 



Figure 4: A plot of the normalization function A\(k) along the equator of a 
sphere in k-space as a function of azimuthal angle, 0. The sharp spikes occur at 
the intersections of the equator with the coordinate axes. At these points the T 7 
hybridization matrix elements vanish exactly, and the width of these peaks is set by 
the ratio of the two Kondo temperatures T o7 /T o8 T o7 j Aqef- The behavior for a 
fixed azimuthal angle (0=0) as a function of the polar angle, 9, is the same. 

Figure 5: Plot of the lowest energy band, E^, for the two mean-field parameter 
sets (labeled as (a) and (b)). In set (a), the conduction electron chemical potential 
is 2.44 eV, and in set (b), it is 3.33 eV. Note that the bands are just shifted with 
respect to each other. For both parameter sets the unshifted T 7 multiplet energy sits 
at -2.0 eV. 

Figure 6: A diagrammatic representation of the infinite-order summation that 
gives the mean-field dressed f-electron propagator. All the Green functions are de- 
fined pictorially at the top of the figure. The bare conduction electron propagator is 
represented by a solid line. 

Figure 7 : Diagrammatic representation of Dyson's equation for the conduction 
Green function. The double line is the hybridization dressed conduction Green func- 
tion, and S co .(/c) is the conduction self-energy. All other elements of the figure are as 
defined in Figure 6. 

Figure 8: Diagrammatic representation of Dyson's equation for the mixing Green 
function. 

Figure 9: (a). Dyson's equation for the diagonal component (in the crystal 
field indices) of the inverse of the (matrix) Bosonic propagator -Dpp. (b). The three 
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hybridization-dressed (mean-field) Fermionic Green functions. Note that, in general, 
the / Green function can mix the multiplet indices at k-points away from the zone 
center. This is because away from the zone center, the symmetry is lower than 
cubic, allowing the crystal field indices to mix. (c). Leading contributions to the 
components (ss, sA, and AA) of the self-energy matrix from closed Fermionic loops. 
The "x" symbol represents the scaled hybridization matrix element, Vr aa (k), which 
is assumed to be of order one. The external legs in the T ss rr diagram come from the 
scaled fluctuations in the s fields, Ss^p. 

Figure 10: Numerically evaluated function I ss (q) as a function of the mesh pa- 
rameter, n mes h for mean-field parameter set (a). Here we used q = 0.5z. 

Figure 11: Numerically evaluated function I s \(q) as a function of the mesh 
parameter, n mes h for mean-field parameter set (a). Here we used q = 0.5z. 

Figure 12: Numerically evaluated function I\\(q) as a function of the mesh 
parameter, n mes h for mean-field parameter set (a). Here we used q = 0.5z. 

Figure 13: Quasiparticle scattering amplitude for incoming particles (solid lines) 

— * — * 

of momenta k and -k. The wavy lines are dressed boson propagators, and the vertices 
denote the anisotropic coupling of quasiparticles to bosons. 

Figure 14: Fermi surface averages in the A lg and T 2g pairing channels of the 
scattering amplitude, (T)/T o7 , as a function of the Boson mesh parameter n mes h- 
& n mesh equals the number of tetrahedra used in the Brillouin zone integrals for the 
Bosonic Green functions. These data are for mean field parameter set (a). (See Table 
5.3.) 

•The horizontal line marked by T 29i / oca z denotes the size of the local ("hard-core") 
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contribution to the average in the T 2g channel. The corresponding local contribution 
for the Ai g channel is 3.54 and would be just above the top of the graph. We see 
that the nonlocal contribution in the T 2g channel is of the same size as the local 
contribution. Even with the fluctuation at n mesh =23, it seems unlikely there is a 
pairing instability of T 2g symmetry. 

•In the Ai g channel, even with the large fluctuation for n mes h=S9, there is clearly no 
pairing instability. 

Figure 15: Fermi surface averages in the E lg and T lg pairing channels of the 
scattering amplitude, (T)/T o7 , as a function of the Boson mesh parameter n mes h- 
These data are for mean field parameter set (a). The solid (dashed) horizontal line 
denotes the local contribution to the average in the E g {T\ g ) channel. 
•Although not yet converged, the averages in the T\ g channel point to the possibility 
of a pairing instability. Averaging the values of {T) T i g /T o7 for all the values of n mes h 
used here yields 

^ Tlg = -0.212 ±0.025. 
T o7 

•In the E g channel, there is a large fluctuation at n mes / l =39. Therefore, we are 
hesitant to say that this is evidence of a superconducting instability. However, it is 
clear that the average interactions in this channel are attractive for a relatively wide 
range of mesh sizes: 25< n mes h < 47. 

Figure 16: Fermi surface averages in the A lg and T 2g pairing channels of the 
scattering amplitude, (T)/T i, as a function of the Boson mesh parameter n mes h- 
(,% n mesh equals the number of tetrahedra used in the Brillouin zone integrals for the 
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Bosonic Green functions.) These data are for mean field parameter set (b). 
•In the T 2g channel, the horizontal, dashed line denotes the local contribution to the 
average. The fluctuations in the full average are at least a factor of two larger than 
the local part and are also varying about zero. Thus it is not possible to say, from 
the present data, if there is a T 2g instability or not. Using the results from parameter 
set (a), however, it still seems unlikely that there is an instability in this channel. 
•In the Ai g channel, it is easy to see that there average interactions are repulsive and 
strong. There is no instability in this channel. 

Figure 17: Fermi surface averages in the E\ g and T\ g pairing channels of the 
scattering amplitude, (T)/T o7 , as a function of the Boson mesh parameter n mes h. 
These data are for mean field parameter set (b). 

•In the Eg channel, the horizontal line denotes the contribution to the average in 
the local limit. For the full average, up to n mes h—Sl, the average value is fluctuating 
evenly about the local value. The surprisingly large fluctuation at n mes / l =31, however, 
makes it impossible to tell if there is an instability in this channel. 
•In the Ti g channel, the full average, although fluctuating, remains negative for 21< 
n mes h < 41. Averaging these values gives a result of 

^ Tlg = -0.226 ±0.049. 
T o7 
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parameter 


V 


So 


£7 




T o7 


D 


set 


(eV) 




(eV) 


(eV) 


(K) 


(eV) 


a 


0.8595 


0.093472 


0.0217648 


0.0208399 


9.25 


2.4405 


b 


0.650 


0.142409 


-0.8437021 


-0.844999 


13.0 


3.33856 
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parameter set 


x 7 


x 8 




Vs 


a 


1.820 


0.0717 


0.975 


0.04029 


b 


1.729 


0.0949 


0.925 


0.0498 



62 



T) 




Alg 


7=2 


Foo + 0.76376261F 4 o + 0.4564355(^4 + I4-4) 


) 


M g 


0.58630197(F 62 + F 6 _ 2 ) - 0.3952847(F 66 + F 6 _ 6 ) 


E g 
E 9 


75(^22 + ^2-2) 


Ti 9 


( -0.93541435i ^ v \r \ 0.353553391i ( v v V 
{ ^2 \1AV I4-1J ^ l y 43 *4-3j, 

( 0.93541435 (\r , \r \ -0.353553391^ 1 V \\ 
{ ^ 1 l4 ~ i ' ) V2 l y 43 + ^4-3jJ 

75(^44 - K4-4) 




T2 9 


73(^21 - Y 2 -i) 
7l(^21+^2-l) 
75(^22 - ^2-2) 
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Representation Dimensionality Transforms like 



Ai 1 x 2 + y 2 + z 2 

E 2 x 2 - y 2 , 3z 2 - r 2 

T\ 3 x, y, z 

T 2 3 xy, yz, zx 
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(a) Cubic Symmetry 


(r \ocal) 'r)/T 7 


^ local) n/Tol 


(r local, hot) /T 7 


V 


set (a) 


set (b) 


set (a) 


A lg (s—wave) 


3.54 


3.53 


-0.215 


Eg(d x 2-y2, d 3z 2_ r 2) 


-0.0275 


-0.0453 


-0.0696 


Tig 


0.0478 


0.0695 


0.0 


T2g{d X y, dyz, d xz ) 


0.0487 


0.0709 


0.0 
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Spherical Symmetry 


(r local) rj/T 


i]=0 (s-wave) 
i]=2 (d-wave) 
77=4 (g-wave) 


1/3 
8/21 
2/7 



